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A dual-scale model of the torrefaction of wood was developed and used to study industrial configura¬ 
tions. At the local scale, the computational code solves the coupled heat and mass transfer and the ther¬ 
mal degradation mechanisms of the wood components. At the global scale, the two-way coupling 
between the boards and the stack channels is treated as an integral component of the process. This model 
is used to investigate the effect of the stack configuration on the heat treatment of the boards. The sim¬ 
ulations highlight that the exothermic reactions occurring in each single board can be accumulated along 
the stack. This phenomenon may result in a dramatic heterogeneity of the process and poses a serious 
risk of thermal runaway, which is often observed in industrial plants. The model is used to explain 
how thermal runaway can be lowered by increasing the airflow velocity, the sticker thickness or by 
gas flow reversal. 

© 2013 Elsevier Ltd. All rights reserved. 


1. Introduction 

In order to prevent attack of wood by fungi or insects, several 
chemical substances are commonly used. However many of these 
are hazardous to the environment. The heat treatment of wood un¬ 
der inert atmosphere at temperature levels ranging from 200 °C to 
300 °C (mild pyrolysis) is an alternative to the use of chemicals, as 
it increases the dimensional stability and durability of wood. Due 
to a thermal degradation of the cell wall constituents, mild pyroly¬ 
sis significantly changes the material properties, which results in 
better durability and stability than natural wood [1-4], However, 
one has to keep in mind that this is to the detriment of decreased 
mechanical strength, density and surface roughness [5,6], In spite 
of these drawbacks, heat-treated solid wood has become popular 
and several production plants now exist throughout Europe [7], 
These industrial plants must overcome serious problems such as 
thermal runaway and internal checking to obtain a relevant and 
homogeneous final product quality. There are several explanations 
for these phenomena, such as the huge variability of biological 
products and the dispersion of the initial moisture content and 
particle shape. However, the dual spatial coupling (within boards 
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and within the stack) in the presence of exothermic reactions is 
likely to be the main source of problems. 

Furthermore, the growing interest in biomass as a source of 
renewable energy gives rise to a crucial demand in biomass pre¬ 
conditioning. For example, torrefaction is a pretreatment required 
in the multistep process from biomass to liquid (BtL), which in¬ 
cludes a gasification followed by a Fisher-Tropsch synthesis [8,9], 
In this case, the reduction of the mechanical resilience becomes 
an advantage, allowing the fibrous behavior to be lost and fine par¬ 
ticles and powder flowability to be obtained with a reduced grind¬ 
ing energy [10]. In this case, the optimal treatment can be tuned 
according to the global mass loss due to thermal degradation 
[11]. However, it is very difficult to obtain this ideal treatment be¬ 
tween and among particles, which emphasizes the need for predic¬ 
tive modeling tools. 

Over the last few decades, numerous works have been devoted 
to the modeling of coupled heat and mass transfer in porous med¬ 
ia: formulation, characterization, computational solution, and val¬ 
idation [12], However, the specific domain of biomass torrefaction 
involves reactive transfer because chemical reactions occur to¬ 
gether with heat and mass transfer. The thermo-activated nature 
of the chemical reactions together with the presence of reaction 
heats explains why this coupling is two-way: the reaction kinetics 
depend on the local value of temperature whereas the temperature 
field is influenced by the sink/source terms due to the reaction 
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Nomenclature 



C p 

specific heat Q kg 1 K ') 

Sc 

Schmidt number 

C 

molar concentration (mol m~ 3 ) 

Sh 

Sherwood number 

Dab 

diffusion coefficient of vapour in air (m 2 s _1 ) 

T 

temperature (K) 

A Y k 

length of the control volume kin the longitudinal direc- 

x v 

molar fraction of vapor 

AZ k 

tion of the board (m) 

height of the control volume kin the thickness direction 

x 

airflow direction (m) 


of the board (m) 

Superscripts and Subscripts: 

h 

external heat transfer coefficient (W m 2 I< 1 ) 



h m 

external mass transfer coefficient (m s ’) 

b 

bound 

j 

mass flux exchanged between the boards and the air¬ 


end board face 


flow (kg s 1 m 2 ) 

j 

general index for a board layer 

Q 

mass flux in channels formed by two adjacent board 

k 

number of the CV on the board face 


layers (kg s 1 m 2 ) 

lo 

lower board face 

L 

characteristic length (m) 

Ncv 

total number of CVs on the board face 

M v 

the molar mass of vapor (kg mol 1 ) 

surf 

surface 

Nu 

Nusselt number 

up 

upper board face 

Pr 

Prandtl number 


vapour 

Re 

Reynolds number 

oo 

external value 


heat. The complete modeling of these phenomena, at the particle 
level or at the level of a particle set, either for a fixed or agitated 
bed, is complex and was addressed here using a stepwise approach. 

Most studies performed using a thermo-balance could be con¬ 
sidered at the micro-particle level, a spatial size small enough for 
the time constant of the kinetic reactions to be larger than all time 
constants related to the heat and mass transfer inside the particles. 
This configuration is easy to obtain at low heating rates, but is 
more delicate to define at very large heating rates, such as those 
encountered in flash pyrolysis [13], 

The literature available at the micro-particle level is much more 
abundant; see for example [14-17], From the literature data, it be¬ 
comes clear that multistage mechanisms are required for the mod¬ 
el parameters to be valid over a certain range of temperature and/ 
or heating rate values. In spite of the numerous papers devoted to 
this subject, one has to emphasize that most of the works, primar¬ 
ily based on the analyses of experimental data, are performed at 
linearly increasing temperature, with heating rates of some de¬ 
grees per minute. This protocol prevents the slow reaction occur¬ 
ring at mild temperature levels to be accurately described. 

At the macro-particle level, the primary variables (temperature, 
moisture content, degradation levels, internal pressure...) evolve 
in time and space. The variable fields must therefore be computed 
as the process evolves. This is the perfect scale at which reactive 
heat and mass transfers are to be accounted for. All models pub¬ 
lished in this field are more or less comprehensive, each one having 
its strengths and weaknesses. For research involving biomass pyro¬ 
lysis, the strong point often lies in the reaction mechanisms and 
the volatiles production. This implies that the overpressure due 
to volatiles production and the induced gas velocities are usually 
taken into account [13,18-21], Some of this research considers 
the change of physical properties and shrinkage from native wood 
to char. However, in most models, less attention is paid to moisture 
migration: either the particle is supposed to be fully dried before 
the process commences, or water evaporation is governed by a 
simplified expression. This is certainly due to the high temperature 
levels encountered during rapid pyrolysis, which causes the parti¬ 
cle to dry very fast. A more recent work [22] proposed a compre¬ 
hensive model able to account for most of the phenomena 
involved at the level of a macro-particle. For example, their model 
perfectly predicts the plateau at the boiling point of water when 
dealing with humid particles. 


However, to our knowledge, none of these works focused their 
attention on the overshoot in temperature due to exothermic reac¬ 
tions. This phenomenon is very important in the domain of mild 
pyrolysis, where the temperature level is low (typically in the 
range 200-300 °C), but the residence time is very long (up to sev¬ 
eral hours). In such conditions, the effect of some tens of degrees 
might have a dramatic impact on the heterogeneity of the heat 
treatment. Based on our extensive knowledge in the field drying, 
we proposed a comprehensive model of heat and mass transfer 
able to deal with any MC field inside the product, including liquid 
water [23], As a result of this feature, our model is perfectly able to 
produce the two pressure peaks measured during the heat treat¬ 
ment of a thick board: the first one being due to water evaporation 
and the second one due to volatiles production. In addition, by 
using both endothermic and exothermic reactions in the selected 
kinetics model, we were able to simulate the temperature over¬ 
shoot observed in experimental data. Using this model, we proved 
that the presence of water delays the heat treatment phase and 
produces a more sudden temperature increase, which exacerbates 
the temperature overshoot. Notice that such an overshoot seems to 
be visible in the experimental data presented in [22] for a 44% 
moist wood log. 

At the upper scale of a set of particles, either fixed or agitated, 
the published modeling work is rare. Di Blasi [24] modeled wood 
gasification in a 1-D fixed bed of particles. This model, which in¬ 
cludes a drying zone and a pyrolysis zone, is able to capture the 
evolution of variables such as the solid and gas temperatures, mo¬ 
lar fractions and reaction rates. Note that, as constant internal 
fields are assumed within each particle, this model is a single scale 
model. 

Compared to this short literature review, the present work pro¬ 
poses a full dual-scale model, whereby each particle is accounted 
for at the upper scale. This gives rise to a full two-way coupling be¬ 
tween the upper model (gas flow) and the particle level (macro¬ 
scopic formulation). The kinetics model is supplied at each node 
of each particle model. Furthermore, the present work uses multi¬ 
stage mechanisms selected from different literature sources for 
each component of lignocellulosic components: cellulose, hemicel- 
luloses and lignins [25]. Even though most literature data were ob¬ 
tained over a large range of temperature levels, the macroscopic 
model was adequately validated for mild pyrolysis [23], The local 
model solves the set of equation in 2D, allowing for the effect of 
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huge anisotropy of wood to be accounted for [26], The dual-scale 
formulation adopted in this work follows the strategy proposed 
for the drying process [27,28], Indeed, for what we believe is the 
first time, we propose a full dual-scale model of the heat treatment 
of wood (a stack of boards or a thick bed of particles), using a com¬ 
prehensive 2-D formulation of reactive heat and mass transfers at 
the particle level, a configuration which is relevant at the industrial 
scale. 

We first provide a brief overview of the dual scale model. 
Among the important features embedded in the comprehensive 
physical and chemical formulation at the particle model, emphasis 
should be given to (i) its ability to deal with any amount of mois¬ 
ture inside the product, including liquid water, (ii) its ability to 
consider water evaporation and production of volatiles on the 
internal pressure and (iii) its ability to predict the temperature 
overshoot due to exothermic reactions. 

At the global scale, the two-way coupling between the boards 
and the channels formed by two adjacent board layers is consid¬ 
ered throughout the process. 

Thereafter, selected simulations are proposed to emphasize the 
effects of the initial moisture content and the board thickness on 
the heat treatment of boards. These effects are studied at both 
the board and stack scales. Thereafter several configurations are 
studied to limit the thermal run-away due to exothermic reactions 
within the board thickness and along the stack. Finally, following 
these test cases, an industrial treatment process is studied. 


2. Dual-scale modelling of wood torrefaction 

2.1. Modeling heat and mass transfer and thermal degradation at the 
local scale (the board) 

At the local scale, the heat and mass transfer and the chemical 
modifications of the cell walls during the torrefaction period are ta¬ 
ken into account in a multiphysics computational model known as 
TransPore. Indeed, during heat treatment, a two-way coupling ex¬ 
ists between the heat and mass transfer (bound water, free water, 
vapor water, volatiles matters) and chemical mechanisms, includ¬ 
ing the thermal activation of the chemical reactions, together with 
the treatment of heat sources (or sinks) and the production of vol¬ 
atiles (Fig. la). The classical macroscopic formulation of heat and 
mass transfer phenomena that arise in hygroscopic porous media 
is presented in detail in [29], The thermal decomposition of cellu¬ 
lose, hemicellulose and lignin components of the wood cell occurs 
via a series of chemical reactions [25]. The reaction rates are as¬ 
sumed to obey an Arrenhius law and each reaction has its own heat 
of reaction, either endothermic or exothermic. The temperature 
computed by TransPore for each control volume and at each dis¬ 
crete time level is used to obtain the reaction rates. Then, the 
resulting heats of reaction and gas production/consumption rate 
are computed by the pyrolysis model and integrated into the heat 
and mass balances computed by TransPore. This comprehensive 
model and its intricate two-way coupling, arising between the cou¬ 
pled heat and mass transfer and the degradation reactions, have 
been validated with experimental data [23], 

Due to the huge permeability ratio between longitudinal and 
radial directions, the board length is of utmost importance in high 
temperature configurations [29], We recently proved that the lon¬ 
gitudinal migration also significantly affects the heat treatment 
[26], Therefore, the 2D version of TransPore was chosen here at 
the particle level. In this work, only half the board length is com¬ 
puted, with symmetry boundary conditions applied on one side 
(see Fig. la). This computational domain is discretized by a set of 
21 x 51 control volumes (CVs) for a 18mm-thick board and 
41 x 51 CVs for a 35 mm-thick board (radial and longitudinal 


directions). In the longitudinal direction, as the computed domain 
is 1 m long, the mesh is refined at the endpiece following a geo¬ 
metrical progression in such a way that the CVs have almost a 
square shape [30], 

2.2. Modeling heat and mass transfer at the global scale (the stack) 


The gas characteristics throughout the stack depend on the 
behavior of each single board and vice versa. A two-way coupling 
arises between the gas flow and the boards during the drying per¬ 
iod [31] and also during the torrefaction period, where the heat 
produced by exothermic reactions occurring in each board can 
propagate over other boards situated downstream. In industrial 
plants, the gas flow is produced by the important pressure gap gen¬ 
erated by fans placed upstream from the stack inlet. The gas flows 
in all tunnels formed between successive lumber layers by stickers 
placed along the flow every 50-70 cm. By construction, the gas 
flow through the stack is therefore simply a collection of mono- 
directional flows acting in parallel. This specific configuration 
allows the evolution of the gas characteristics (drying and heat 
treatment conditions throughout the board stack) to be captured 
using a simplified formulation. The energy and mass balances are 
written in each channels formed by two adjacent board layers 
(refer to Fig. lb and c). In the case of heat treatment, volatiles pro¬ 
duced by the chemical reactions are assumed to be perfectly mixed 
with the dry gas (N 2 ) of the gaseous phase assumed to form a 
homogeneous gas mixture [23], In this way, the mass flux of this 
dry mixture, j a exchanged between the boards and the airflow 
has to be taken into account in the heat and mass balances (see 
Eqs. (2) and (3)). C pa in Eq. (3) is the specific heat at constant pres¬ 
sure of the dry mixture, which is assumed to be equal to the value 
of dry gas. In the following equations, each airflow channel is sur¬ 
rounded by upper board faces (superscript up of board layer j) and 
by lower board faces (superscript lo of board layer j + 1). 

Moisture balance 


_ ;upj , Joj+l 

dx +lv 
Balance of dry mixture 
dq a , upi ., 0J+1 
dx =J ° +} ° 

Enthalpy balance 

(C pa q a + C pv q v .) ^ = T +ic m +jy J (C pv (TZ - T) 

+f v j+ \Cnv(T , $' ~ T) +ja J+ '(C pa (T , $' 

+j7 J (C pa Cl%i-T) 


( 1 ) 

(2) 


T) 

( 3 ) 


Note that the binary diffusion in the airflow is neglected due to the 
high value of the Peclet number. In the conservation Eqs. (1)—(3), 
the heat flux (j c ), the mass flux of vapor water (j v ) and the mass flux 
of dry mixture (j a ) are provided by the solution of each particle (via 
computation using one module of TransPore for each particle). In re¬ 
turn, the solution of Eqs. (1)—(3), allows the partial vapor pressure 
and temperature of the airflow to be computed throughout the 
whole stack and to be used thereafter as boundary conditions in 
each module of TransPore. In this work, the mass flux of water vapor 
and dry mixture coming out of each board face are obtained as a 
weighted average over the CV faces of that board, including half 
the flux at the end piece. For example, the flux at the upper face 
of the board layer j is obtained using the following formula: 


MP J7k A yk + E*i~ end -x- tek 

ESTA y k 


(4) 
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(b) 



~>r~ 

A stack of boards 


Boards layer 


(C) 


Local fluxes leaving 
each face of the 


R U 



The usual boundary layer theory is used to estimate the external 
heat and mass transfer between the surface of the boards and the 
gas flow. The fluxes of total moisture and total enthalpy at the 
boundary respectively are of the following form: 

L = h m cM v In (\——^ 

' (5) 

Jc = HTb - T„) 

The gas flow channel between two adjacent board layers can be 
viewed as flow between parallel plates. The heat transfer coefficient 
is calculated from the Dittus-Boelter equation [32]: 

Nu = 0.0230Re 4/5 Pr 1/3 (6) 

where Nu is the Nusselt number, Re the Reynolds number and Pr the 
Prandtl number. The characteristic length L used to estimate the 


Nusselt and Reynolds numbers is equal to the hydraulic diameter, 
defined as four times the cross section area divided by the wetted 
perimeter, and taken as twice the sticker thickness in our 
configuration. 

Eq. (6) is considered valid for turbulent flow (Re > 10,000). In¬ 
deed, the gas flow inside an industrial plant is highly turbulent, 
due to the presence of fans, the surface roughness of sawn boards 
and the geometrical defaults between successive boards. Conse¬ 
quently, Salin [33] suggested that Eq. (6) can be used with accept¬ 
able accuracy for lower values of Re. This turbulent regime also 
explains why constant transfer coefficients can be used throughout 
the stack. 

The rectangular duct between two board layers and two stickers 
has a large dimensional ratio (more than 50). Thus, in this config¬ 
uration, the hydraulic diameter can be estimated at twice the 
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External heat and mass transfer coefficients calculated according to the airflow 
velocity for 2 cm thick stickers. 


Gas velocity (m/s) 2 5 8 11 

h (W nr 2 °C -1 ) 14 29 42 55 

ft m (ms-') 0.014 0.028 0.041 0.054 


sticker thickness. The airflow velocity and the external heat and 
mass transfer coefficients are assumed constant along the air flow 
path and across the duct section. The external mass transfer coef¬ 
ficient is calculated using the analogy between heat and mass 
transfer: 

h m = D M Sh/L (7) 


code TransPore for the reference board (A). 


Thickness x Length 


Basal density 
Cellulose content 
Hemicellulose content 

Gaseous permeability in radial direction 
Gaseous permeability in longitudinal direction 
Diffusion resistance in radial direction 
Diffusion resistance in longitudinal direction 
External mass transfer coefficient 
External heat transfer coefficient 


Sh = Nu■ Sc ,/3 Pr 1/3 


Table 1 depicts the estimated external transfer coefficient for a 2 cm 
thick sticker and according to the different airflow velocity along 
the stack. Additionally, Table 2 depicts the estimated external trans¬ 
fer coefficient for a 2 m/s airflow velocity for different sticker 
thicknesses. 


Table 4 

Input values for the reference stack of boards used at the global scale. 

Value 

Board dimensions (width x length) 200 x 2000 

Number of boards 11 

Sticker thickness 2 


2.3. The test configuration 

As product parameters such as thickness and initial moisture 
content dramatically affect the homogeneity of the final product, 
3 types of beech boards are simulated in the following section: 

• Case-study A: a reference particle (3% MC, 35 mm thick) with 
parameter values as specified in Table 3. 

• Case-study B: Higher moisture content (20% MC, 35 mm thick). 

• Case-study C: Higher moisture content with a thinner board 
(20% MC, 18 mm thick). 

For each case study, a full layer of boards has been computed to 
analyze the effect of the stack. The stack is composed of 11 boards 
arranged in layers (thus, 2.2 m deep along the air flow). The stack 
configuration is summarized in Table 4. In order to focus our atten¬ 
tion on the stack effect rather than the heterogeneity of the initial 
product, all simulations were performed with identical boards in 
the present work. This allowed us to limit the stack calculations 
to one single stack layer. 

The following schedule is applied as the heat-treatment 
procedure: 

• Heating up to 110 °C at a rate of 2 °C/min, followed by a 5-h pla- 
teau (drying phase), 

• Heating up to 220 °C at a rate of 1.5 °C/min, followed by a 5-h 
plateau 

• Cooling down to 90 °C at a rate of -1.5 °C/min. 

Using a personal computer with a 2.93 GHz i7 Intel processor, 
the CPU time is about 3 min for one single board simulation 
(Case-study A) and rises to about 40 min for a full board layer 
simulation. 


Table 2 

External heat and mass transfer coefficients calculated at 2 m s~ for different sticker 
thicknesses. 


Sticker thickness (cm) 2 5 8 11 

h (Wm -2 ”^ 1 ) 14 12 11 10 

h m (rns 1 ) 0.014 0.011 0.010 0.010 


3. Results and discussion 

3.1. Coupled mechanisms arising due to the dual-scale (case study A) 

Our first comment will be focused on case study A, the thick 
board (35 mm) with almost no initial moisture content (3% MC). 
Fig. 2 depicts the evolution of the temperature at surface and core 
(half-length) and the dimensionless pressure at the core for boards 
1, 6 and 11, respectively situated at the stack inlet, in the middle 
and at the stack outlet. 

Due to the one-way interaction along the gas flow, the first 
board is simply submitted to the process conditions set by the 
above-defined schedule and is not influenced by the presence of 
the stack. At first, we comment on the results obtained for this inlet 
board to provide a good understanding of the coupled mechanisms 
arising at the particle level, before advancing further towards the 
stack effect. A classical convection/conduction temperature profile 
is observed during the first temperature increase (up to the drying 
plateau at 110 °C). The surface temperature increases at first, as a 
result of the heat supplied by the gas flow by convective exchange 
and the core temperature follows this rise with a time delay due to 
conductive heat transfer. This transient period lasts for around 2 h 
and the temperature field is almost perfectly equilibrated at 110 °C 
after 3 h. 

A second transient period appears at the end of the drying pla¬ 
teau, when the gas temperature rises up to the heat treatment pla¬ 
teau. As before, a temperature difference develops between the 
gas, surface and core temperatures, as a result of the resistances 
to heat transfer placed in series (convective transfer in the bound¬ 
ary layer and conductive transfer inside the product). Then, when 
the temperature overpasses 200 °C, one can see the core tempera¬ 
ture increasing more quickly than the surface temperature, and 
eventually overtaking it. This is obviously the effect of exothermic 
reactions, triggered by the temperature level, which produces an 
internal source term. A temperature overshoot is observed, during 
which the convection/conduction mechanism act in the opposite 
way, the core being warmer than the surface and the external heat 
exchange coefficient operating to cool down the product. A maxi¬ 
mum core temperature of 240 °C is observed for the inlet board 
after 3 h at the 220 °C-plateau. 
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The double pressure peak is always very informative regarding 
the sequential occurrence and intensity of drying and heat-treat¬ 
ment. The first gas overpressure is due to water vaporization inside 
the product when the temperature significantly increases above 
the boiling point of water and the second peak appears due to 
the volatiles produced by the decomposition reactions [23], In case 
study A, the first peak is weak because the initial MC is low. The 
second pressure peak reaches 1.7 bar for the first board. 

Generally speaking, when the gas flow has to supply heat to the 
board, the gas flow temperature decreases along the stack and, 
progressively, the gas flow losses its ability to supply heat. On 
the contrary, when the boards supply energy to the airflow due 
to the exothermic reactions, the gas temperature increases along 
the stack and its ability to cool down the product worsens. In addi¬ 
tion, the increased gas temperature is able to further activate the 



Fig. 3. Evolution of the gas flow temperature throughout the stack (case study A). 


degradation reactions. The same trend holds for moisture ex¬ 
change. For example, during the drying phase, the dew point in¬ 
creases along the stack due to the vapor flux exiting the boards. 
Together with the temperature drop, this effect further decreases 
the drying capacity of the gas flow. To summarize, the stack effect 
amplifies the trends observed for one single board. The interest of 
our dual-scale model lies in its potential to quantify these accumu¬ 
lative effects. 

The initial heating up period perfectly confirms the above-men¬ 
tioned trend: the time delay between the gas flow temperature 
evolution and the surface and the core temperatures of the board 
increases along the stack: 4 h (increasing period + plateau) is re¬ 
quired for board 11 to approach the dry bulb temperature, com¬ 
pared to 2 h for the inlet board. 

The evolution in time of the gas temperature profile along the 
stack (Fig. 3) perfectly depicts the existence of phases during which 
the gas temperature decreases from inlet to outlet (typically the 
two heating up periods) and phases with the opposite effect (the 
period of intense chemical degradation and the cooling down per¬ 
iod). At the end of the drying plateau, the gas temperature profile is 
flat, which denotes a period without any exchange between the gas 
flow and the boards (one has to recall that Fig. 3 is plotted for case 
study A, which has almost no moisture to be removed). 

The delay in the heat treatment of the outlet board (n°ll) in 
comparison with the inlet board can be assessed by the second 
pressure peak, which occurs one hour later than for board 1 
(Fig. 2). The higher level of temperature at the stack outlet intensi¬ 
fies the magnitude of the temperature overshoot. A thermal run¬ 
away is evident for board 11, whose temperature attains 266 °C 
at its core and 252 °C at its surface and with a pressure peak at 
2.23 bar. A detailed analysis of the temperature field inside the 
board tells us that the exothermic reactions are first triggered close 
to the endpiece, which heats up first due to the presence of orthog¬ 
onal exchange faces (Fig. 4). The temperature overshoot thereafter 
spreads throughout the length, whatever the board position. Note 
that at the end of the process, when the dry bulb temperature of 
gas at the stack inlet is 90 °C, the internal temperature of the outlet 
board is still at 140 °C. The pressure curves of boards 1, 6 and 11 
confirm along the stack: a slower drying, with wider and lower 
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Fig. 4. Carpet plot of the temperature field inside board 1 (left) and board 11 (right), when the exothermic reactions start and depict a higher value near the endpiece (top) 
and when the core temperature reaches its maximum (bottom). 


vaporization peaks and a treatment period becoming more intense, 
with sharper and higher volatiles peaks. 

3.2. Effect of initial MC and thickness (case studies B and C) 

Case study 2 (20% MC, 35 mm thick) allows the effect of initial 
moisture content to be analyzed (Fig. 5). Due to the high initial MC, 
the drying phase is far from complete at the end of the first plateau. 
This can be observed by the residual temperature between gas and 
surface at 360 min, even for board 1. This is required to supply, by 
convective transfer, the latent heat of evaporation needed to re¬ 
move the moisture as water vapor. Consequently, a large part of 
evaporation takes place during the second heating-up period and 
even at the beginning of the heat treatment plateau. This delays 
the increase of the surface temperature, just because the convec¬ 
tive heat exchange is used to supply the latent heat of evaporation, 
but delays even more the increase of the core temperature, which 
dries last. At this temperature level, drying occurs by internal 
vaporization which gives rise to an overpressure inside the product 
and a longitudinal water vapor migration in Darcy’s regime due to 
the huge anisotropy ratio [29]. The massive and sudden evapora¬ 
tion occurring for board 1 during the second heating-up period, 
and even during the beginning of the heat treatment plateau, gives 
rise to an important pressure peak (around 1.8 bar).An unexpected 
consequence of this drying behavior is that, once the drying of 
board 1 is completed, the temperature gradient inside the board 


is huge. Consequently, the core temperature increases rapidly to¬ 
wards the treatment temperature and triggers the exothermic 
reactions. Interestingly, this board, which requires a much higher 
amount of energy than case study A before the heat treatment pre¬ 
sents a similar value of the maximum core temperature (240 °C). 
This maximum is however attained at 630 min against 580 min 
for case study A. 

During the drying period, the decrease of dry-bulb temperature 
along the gas flow is larger than for case study A (22 °C at 120 min 
against 7 °C previously).The softer drying conditions imposed to 
boards 6 and 11 induce a delay in the first gas pressure peak in 
comparison with board 1 (by about 65 min for board 11). Fig. 6 
proves that board 1 still has a large amount of bound water at 
the beginning of the treatment plateau and is fully dried after 
one hour at this temperature level (480 min). With the stack effect, 
the MC extraction is even longer for the outlet board (n°l 1). At this 
position, drying hardly started at the end of the “drying” plateau 
and a significant amount of bound water still exists after one hour 
at the treatment plateau. Note the wider and lower pressure peaks 
obtained when the drying conditions softens, which is even more 
obvious than for case study A. Such conditions also induce a delay 
in the heat treatment of the outlet board in comparison with the 
inlet board, which can be assessed by the second pressure peak 
occurring 105 min later for board 11 with maximum temperature 
levels of 264 °C at the core and 244 °C at the surface. In this case, 
the temperature overshoot seems to be lower than for the case 
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Fig. 5. Evolution of temperature (surface and core) and dimensionless pressure (core) for boards 1, 6 and 11 (from inlet to outlet). Case study B (MC = 20%, 
Thickness = 35 mm). 


study A, but this is only because the time delay is so important that 
the treatment occurs when the gas flow already started to de¬ 
crease, hence cooling down the board from its surface. 

In the case of a thin board (case study C, Fig. 7), the time con¬ 
stants inside each board are reduced; typically by a factor 4 for 
all diffusive phenomena. In addition, as the water mass and the 
heat capacity are proportional to the thickness, the stack effect is 
typically reduced by a factor 2.This simple analysis is able to ex¬ 
plain all the trends observed with 18-mm thick boards compared 
to 35-mm thick boards (case study B): 

• The difference between the surface and core is now hardly vis¬ 
ible, at least on the temperature values, 

• The time between the two pressure peaks is significantly 
reduced; 

• The cumulative effect of exothermic reactions is still significant, 
with a maximum temperature increasing from 227 °C for board 
1 to 246 °C for board 11 and an overpressure peak from 1.28 to 
1.50 bar. 

Therefore, case study C tells us that, among all coupled mecha¬ 
nisms observed during this dual-scale configuration of wood torre- 
faction, the most difficult effect to reduce is the thermal run-away 
of the exothermic reactions due to the stack effect. This can be ex¬ 
plained by the thermal activation of these reactions, which are 
therefore significantly faster as soon as the gas temperature in¬ 
creases due to the stack effect. In turn, this activation further in¬ 
creases the gas temperature for the next boards. This coupling is 
therefore likely to appear even for small particles, like those 
encountered in the use of lignocellulosic materials for energy 
purposes. 

Fig. 8 depicts the density field of cellulose at the end of the ther¬ 
mal treatment process for different configurations, at the stack in¬ 
let and outlet. These graphs allow both the heterogeneity inside 
the particle and among particles to be quantified. In our example, 
the worse case for the outlet board is by far the reference configu¬ 
ration (A). The highest level of thermal degradation arose at the 
core of the outlet board, due to the cumulated effects of the mac¬ 
roscopic and global scales. Surprisingly, the treatment seems to 


be more homogeneous for the wet board (B): we just have to keep 
in mind that this effect is solely due to the important delay: for the 
outlet board, the heat treatment occurred at the end of the sche¬ 
dule and the product was already cooled down by the gas flow. 
This observation opens interesting possibilities for a control strat¬ 
egy able to produce a quite homogeneous treatment despite the 
large thickness. Finally, reducing the board thickness is very effi¬ 
cient in obtaining a homogeneous treatment, including within 
the worst board, the one situated at the stack outlet. 

In summary, the behavior of each board of the stack results 
from coupled reactive transfer occurring inside each particle and 
between particles. 

At the local scale, the heterogeneity of heat treatment observed 
within a single board results mainly from: 

• the overheating of the board due to the exothermic reactions 
inside the product together with the insulation role of wood. 
The effect is significant for thick particles, 

• the heterogeneity in initial moisture content. A higher initial 
MC delays the heat treatment that may therefore occur at dif¬ 
ferent gas temperatures for a multistep treatment protocol. 

At the global scale, the coupling between different boards along 
the gas flow amplifies the effect of heating or cooling of the prod¬ 
uct. Therefore, the transition from the first status (the gas flow 
heats up the product) to the second (the gas flow cools down the 
product) is more sudden in time and huge in amplitude. Degrada¬ 
tion reactions are therefore triggered very rapidly and more inten¬ 
sely, as the gas temperature at the stack outlet can be higher than 
the set value (by almost 30 °C for our case study A). All together, 
this explains why the outlet board shifts very rapidly from a slow 
drying/heating phase to an intense treatment phase. 


4. Reducing the cumulative effect of exothermic reactions 

In order to exemplify the practical potential of this sophisti¬ 
cated dual-scale model, this section summarizes the effect of some 
process conditions on the computed temperature overshoots. 
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Fig. 6. Carpet plot of the moisture content field inside board 1 (left) and board 11 (right), at the end of the "drying” period (345 min), at the beginning of the treatment 
plateau (420 min) and after one hour of heat treatment (480 min). 


Different strategies are studied for reducing the treatment 
heterogeneity: 

• increasing the airflow velocity, 

• increasing the sticker thickness, 

• applying different periods of airflow reversal. 


4.3. Increasing the airflow velocity 

Increasing the airflow velocity is a standard approach taken to 
reduce the stack effect, by decreasing the variation in the variable 
fields along the stack. This must have a positive effect regarding 
the treatment homogeneity. On the other hand, a large airflow 
velocity increases the transfer coefficients (see Table 1), which 
tends to deepen the moisture and temperature profiles inside each 
board. A first observation is that this secondary effect works 
against the treatment homogeneity inside each particle. 


In order to confirm these trends, three additional airflow veloc¬ 
ities were computed in addition to the reference case study (A). 
Fig. 9 shows that by increasing the airflow velocity, the heteroge¬ 
neity of treatment within and between the boards decreases. In¬ 
deed, the most important effect of increasing the gas velocity is 
to reduce the evolution of its properties along the stack. 

Additionally, a higher velocity generates a lower external heat 
and mass resistance, which in turn increases the efficiency of the 
heat exchange between the airflow and the boards. The boards 
can be more easily cooled down by the exchange surfaces during 
the exothermic period and the thermal runaway can be avoided 
to a large extent. The difference between the maximum core and 
surface temperatures of each board is plotted in Fig. 9b. This infor¬ 
mation is interesting because it is strongly tied to the treatment 
heterogeneity inside that board. Its interpretation is subtler, as 
the value depicted for the inlet board seems to be independent 
of the airflow. Indeed, as the velocity increases, the transfer coeffi¬ 
cients at the surface also increase, which would produce a deeper 
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Fig. 10. Effect of sticker thickness: (a) Maximal core temperature reached during the process by the 1st, 6th and 11th board from the inlet side for different thicknesses, (b) 
difference of maximal temperature between the core and surface of the same boards. 


temperature profile inside the board in the absence of an internal 
source term. However, during the exothermic period, the internal 
temperature profile is established to evacuate the heat source. If 
this heat source is almost the same, the temperature difference re¬ 
mains identical. In this case, a larger external heat transfer coeffi¬ 
cient will simply reduce the temperature values by a constant 
amount. The maximum temperature will therefore be reduced at 
larger velocities, which is confirmed by Fig. 9a. 


4.2. Increasing the sticker thickness 

Increasing the sticker thickness is another solution to decrease 
the stack effect. Here, we increased the sticker thickness by the 
same amount as the gas flow mentioned in the previous paragraph. 
However, by reducing the stack effect, we now also reduce the 
transfer coefficients between the gas flow and each board (see Ta¬ 
ble 2). Lower external heat transfer coefficients reduce the thermal 
coupling between the gas flow and the boards, therefore reducing 
the potential to cool down the exchange surfaces during the exo¬ 
thermic period. As a result, increasing the sticker thickness is less 
efficient at increasing the velocity to limit the thermal run-away 
(Fig. 10a, to be compared to Fig. 9a). In spite of this loss of effi¬ 
ciency, lower transfer coefficients tend to decrease the tempera¬ 
ture gradients inside each board. The temperature difference 
between the core and surface is then slightly lower in comparison 


with the strategy that increases the airflow velocity (Fig. 10b, to be 
compared to Fig. 9b). 


4.3. Airflow reversal with different period of inversion 

In conventional drying it is well known that the airflow reversal 
is a simple means to reduce the stack effect. This principle is there¬ 
fore likely to apply for the heat treatment. In this section, we com¬ 
puted a series of simulations in which the gas flow is inverted at 
different periods of time (no reversal, every hour, every 25 min, 
every 5 min, and every minute). Fig. 11 summarizes the computed 
results in terms of the maximum temperature attained for three 
chosen boards (1 = inlet, 6 = middle and 11 = outlet) and the max¬ 
imum difference between their core and surface temperatures. 

As expected, the gas flow reversal systematically has a positive 
effect on the outlet board, the worst situation regarding the stack 
effect. However, the price to pay for this improvement is an in¬ 
crease of the treatment heterogeneity for the inlet board. The mid¬ 
dle board seems to always slightly benefit from gas flow reversal. 
The effect of the reversal period is not so obvious to analyze. In¬ 
deed, for very short periods of time, the averaging effect over the 
stack is uniform: when the gas flow is inverted every 5 min, inlet 
and outlet boards become similar as their situation tends to be per¬ 
fectly symmetrical. The middle board is the best in terms of tem¬ 
perature overshoot, as it never undergoes the full accumulative 
effect of exothermic reactions along the stack. For larger reversal 



ConfA Period^ Period= Period= Period^ 

(Reference) 60min 25min 5min lmin 



Fig. 11. Effect of gas flow reversal: (a) maximal core temperature reached during the process by the 1st, 6th and 11th board from the inlet side for different reversal periods, 
(b) difference of maximal temperature between the core and surface of the same boards. 
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periods, the symmetry is broken, as the history of each single board 
strongly depends on the sense of the gas flow when the risk of 
thermal run-away is maximal. This explains for example why a 
one-hour period is not beneficial for the inlet board. In this case, 
the period of intense exothermic reaction throughout the stack oc¬ 
curred while the gas flow was inverted. 

In conclusion, the dual-scale model confirms that the gas flow 
reversal is a good way to reduce the heterogeneity due to the stack 
effect, and tells us that, in order to efficiently reduce the treatment 
heterogeneity, a short period of airflow reversal has to be chosen, 
typically 10 min or less. 

5. Conclusion 

A dual-scale computational model of the drying and heat treat¬ 
ment of wood was developed and used to study industrial config¬ 
urations. At the local scale, the computational code solves the 
coupled heat and mass transfer and the thermal degradation 
mechanisms of the wood components. Among the important fea¬ 
tures of this particle model, emphasis should be given to (i) its abil¬ 
ity to deal with any amount of moisture inside the product, 
including liquid water, (ii) its ability to consider water evaporation 
and production of volatiles on the internal pressure and (iii) its 
ability to predict the temperature overshoot due to exothermic 
reactions. At the global scale, the two-way coupling between the 
boards and the channels formed by two adjacent board layers is 
considered throughout the process. 

The exothermic reactions occurring in each single board can 
propagate along the stack and pose a serious risk of thermal run¬ 
away. In our reference configuration, the maximum temperature 
reached at the core of the board increases from 239 °C (inlet 
board), to 264 °C (outlet board). A combination of a high gas flow 
velocity together with an inversion of the flow over a short period 
of time is a good way to dramatically reduce the stack effect. 

We believe that based on its predictive ability, this sophisti¬ 
cated model has great potential, either for improving existing situ¬ 
ations or in testing innovative process control, such as using a non¬ 
constant treatment temperature, which is the most likely way to 
reduce the heterogeneity in treating thick boards. We also intend 
to apply this model to the torrefaction process of a thick bed of par¬ 
ticles, a configuration of particular interest for the use of lignocel- 
lulosic products as source of energy. 
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